Selective-pivot sampling of radial distribution functions in asymmetric liquid mixtures 
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; 

We present a Monte Carlo algorithm for selectively sampling radial distribution functions and 
, effective interaction potentials in asymmetric liquid mixtures. We demonstrate its efficiency for 

■ hard-sphere mixtures, and for model systems with more general interactions, and compare our 

simulations with several analytical approximations. For interaction potentials containing a hard- 
sphere contribution, the algorithm yields the contact value of the radial distribution function. 



I. INTRODUCTION 



Liquid mixtures have been studied in the past decades with objectives ranging from understanding basic theoretical 
properties to answering questions of technological relevance. The simplest model of liquid mixtures, binary hard 
. spheres, has been an important test bed for experimental, analytical and numerical techniques. 
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FIG. f: Asymmetric hard-sphere mixture (left), and its effective fluid of solutes (pair interaction Vy and triplet interaction 
■ Vkim shown, right). 




• ■ In asymmetric binary mixtures consisting of few large particles (solutes) and many small particles (solvent) (see 
Fig- UK one may in principle integrate out all the degrees of freedom of the solvent, and arrive at an equivalent effective 
fluid of only the l ar ge particles, with complicated effective interactions involving pairs, triplets, and, more generally, 
n-tuples of solutes Usually, all contributions beyond the pair interaction are neglected (see, for example, [2fl), 
even though they are quite important [H, 0, [![. Moreover, integral-equation theories or density functional theories [6j, 
among others, allow to study the mixtures directly. For asymmetric hard-sphere mixtures, various approaches based 
on thermodynamic criteria yield empirical equations of state (see for example Q). However, many theories may fail 
in extended regions of parameter space because of the different sizes of solutes and solvent (for problems with the 
integral-equation theories, see P). 

Numerical simulation is a crucial tool to validate the above methods. However it also meets difficulties for asym- 
metric mixtures. On the one hand, the displacement of each solute is highly constrained by the surrounding solvent 
(see @). Pivot-cluster Monte Carlo algorithms [l(| have partially overcome this problem for hard-sphere models 
and for interacting mixtures [12, G3 > as l° n g as * ne overall density of the system is not too high [b| . On the other 
hand, the difference in size of solutes and of solvent particles generates a sampling problem for observables such as the 
pair-distribution functions, which vary strongly on the very small scale of the solvent. To describe these observables, 
one needs many data (generally contained in a histogram with a fine grid). On close approach of two solutes in 
asymmetric mixtures, the pair correlation functions rise steeply. Therefore, the contact value of the pair-correlation 
function cannot be obtained precisely by extrapolation [l5[ . 

In this paper, we introduce a sclcctivc-pivot sampling algorithm, which allows us to compute the ratio of the radial 
distribution functions g(R)/ g(R') for two arbitrary inter-particle distances R and R' . For hard-sphere mixtures, we 
compute the pair-correlation function of the solutes and especially its contact value. We also determine the potential 
of the mean force. We furthermore discuss the extension of the method to systems with arbitrary pair interactions. 
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II. 



SELECTIVE-PIVOT SAMPLING 
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FIG. 2: Selective-pivot move 1 — » 1' (configuration i, left) and return move 1' — > 1 (/, right). In i, the distance R12 = Ri 
Ri'2 = Rf- All pivots resulting in a new distance Rf lie on a sphere Si with center Oi and radius Rf/2 (crosses). 



and 



For simplicity, we first discuss selective-pivot sampling for a system of an isolated pair of three-dimensional hard 
spheres. To sample the ratio of the pair correlation function for two arbitrarily distances, Ri and Rf, we consider a 
configuration i, with a distance Ri between the two particles, and all moves of one particle (1 — > 1') such that the 
new distance is Rf (see Fig. [2]). Any such move can be interpreted as a reflection around a pivot located on a sphere 
Si with radius Rf/2 and center Oi, half-way between the two particles. We choose the pivot randomly on the sphere 
S] (using Gaussian random numbers, see !). The a priori probabilities for proposing the move i — > /, and the return 
move / — > i are, respectively, 

/) = 5- « i . -*(/-► = — OC — S . 

Si Rf Sf H l 

The detailed-balance condition connects the stationary probabilities ir(i) and 7r(/) with the a priori probabilities: 

7r(l)A(i -H- /) = 7r(/).A(/ -H. i). 

It follows that the probabilities of the two configurations are biased by a geometric factor with respect to the constant 
hard-sphere probabilities, 

Si 



For spherically symmetric systems, the pair correlation function depends solely on the distance between particles and 
reduces to the radial distribution function g(R), which is linked with the probability n(R) to observe a particle in an 
infinitesimal spherical shell of radius R, the other one being at its center, by 

tt(R) oc R 2 g{R)dR. 

The biasing factor in cq. flj) cancels the phase-space factor of the radial distribution function, and the ratio of the 
probabilities and 7r(/) equals the ratio of the radial distribution functions: 

5/' 

The above algorithm and the relation of eq. ^ between the probabilities of observing the distances Ri and Rf remain 
valid for two or more solute particles in a box with periodic boundary conditions, in the presence of other components 
(solvent). It suffices to permanently tag two solutes. One Monte Carlo move involves a single pivot, but may transform 
many particles (see |9(), including both tagged ones, so that the distance between them may not change. 
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III. APPLICATIONS 



To validate the selective-pivot sampling algorithm, and to explore possible applications, we compute the distribution 
function of the solutes for asymmetric binary mixtures, and compare it with simulation data obtained with standard 
cluster simulation methods, as well as with analytic approximations. In this system, we also compute the contact value 
of the distribution function. This can be done directly, without extrapolation. Furthermore, we use the selective-pivot 
sampling algorithm to compute the potential of the mean force, that is, the distribution function for two solutes in a 
large bath of solvent, approaching the regime of infinite dilution. 
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A. Radial distribution function of solutes 



We consider an asymmetric binary mixture of spheres with diameters D s (small particles — solvent) and D\ (large 
particles-solute), and size ratio R = D\/D s = 10, at packing fractions rj s = 0.126 and r)i = 0.121, where the cluster 
algorithm performs well (the packing fraction of the solvent, r/ s , is linked to the number of solvents N s in the simulation 
box of volume V by rj s = f(N s /V)D^, etc.). 




FIG. 3: Radial distribution function of the solutes for hard spheres with R — 10, r/s = 0.126 and rj\ — 0.121. Symbols: Cluster 
algorithm (filled squares: standard sampling method; empty squares: selective-pivot sampling); full line: density functional 
theory; dashed line: RHNC-MSA || 



In Fig. [3l we compare the radial distribution function of the large spheres gn(r) obtained from selective-pivot 
sampling and from a standard cluster algorithm and find very good agreement. The contact region and, especially, 
the contact value g\\(r = D\) become accessible only with our sampling scheme. With standard grid sampling, they 
can be obtained only by extrapolation because the probability of observing two solutes exactly at contact is zero. 
Our sampling approach gives sensibly more accurate contact values than the standard method (see [IB]). I n Fig- El 
9\\{D\) = 13.03 ± 0.26. As a common use of simulation data, we also compare the pair-correlation function with 
density functional theory, which agrees very well [UGH) an d with the Ornstein-Zernike integral equation with the 
RHNC-MSA closure, which overshoots in the contact region Q. 



B. Potential of the mean force 



As another application of selective-pivot sampling, we consider the potential of the mean force, the solute pair- 
correlation function in the limit of infinite dilution of the solute particles: 

lim Qn(r,V h ) = exp [-/3CT(r)] , (3) 

where U(r) = u\\(r) + $ ff (r) is the total pair interaction between the solutes, un(r) the direct interaction, & e fi( r ) 
the pair potential of the mean force and r\h the bulk solvent packing fraction. The selective-pivot sampling algorithm 
allows us therefore to compute the potential of the mean force of the effective fluid obtained by integrating out the 
solvent particles of the binary mixture. One simply puts the two tagged solutes into a bath of solvent and determines 
the potential differences for several pairs Ri, Rf. The size of the bath plays almost no role in this simulation because 
the two large particles are kept at two fixed distances, and cannot escape to infinite separation. The simulation can 
thus take place in a large box, approaching the infinite dilution of the two solute particles. 

To illustrate this point, we computed the potential of the mean force between two hard-sphere solutes in a bath 
of hard-sphere solvent with size ration R = 5 and bulk packing fractions rjh = 0.1 and rjt, = 0.2. We considered here 
N s = 5280 and N s — 10560 respectively in a cubic box of L = 30D S . At a difference of other simulation methods 
for computing the effective interaction, this method needs no extrapolation [l9| and has no adjustable parameter 
[20| . The Monte Carlo data agree very well with the potential of the mean force obtained from the Ornstein-Zernike 
equation with the Rosenfeld fundamental measure closure [2l|, [22j (see Fig. |4]). 
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IV. SELECTIVE-PIVOT SAMPLING FOR MODELS WITH GENERAL PAIR INTERACTIONS 

Our sampling method can be directly generalized to more general interactions, which can also be handled by the 
pivot-cluster algorithm (see [13, As an illustration we study the effect of solvation forces in a binary mixture of 
colloids with size ratio R = 5, and a Yukawa tail for unlike pairs of colloids added to the hard-sphere interaction: 

M _ / oo r < As fA s 

Uls[r) ~ \ -eh cxp{-z ls (r - D ls )}/r r > D ls ' W 

m 
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FIG. 5: Pair distribution function of solute particles for R = 5, rj s = 0.11 and r/\ = 0.12. Empty symbols: standard sampling 
method; filled symbols: selective-pivot sampling; squares: hard-sphere; circles: Yukawa interaction with z* s — 10; triangles: 
with z* s — 2.5. 

We used a value of O.^ksT for the depth of the attractive well, which correspond to e* s = e\ s D s / (k^T) — 3/2, and 
considered two values for the inverse range of the attraction force, z* s — z\ s D s — 10 and z* s = 2.5. In Fig. [5] we see that 
data obtained from selective-pivot sampling and from standard cluster algorithm agree very well. The advantages 
and the applications mentioned in the previous section carry over. Fig. [5] shows that the contact value of the radial 
distribution function decreases when the range of the attraction between small and large particles increases, because 
the solvation forces create a thin shell of small colloids which counterbalances the depletion between the large ones. 
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V. CONCLUSION 



In this paper, we have presented a selective sampling cluster algorithm which allows to obtain the radial distribution 
function, including the contact value, and also yields the pair potential of the mean force between two solutes due to 
the presence of the solvent. The algorithm solves a sampling problem for observables with very rich structure on a 
small scale, but does not overcome the remaining limitation of cluster algorithms for liquid simulations, namely the 
restriction to moderate densities. It will be interesting to see whether this restriction can also be overcome in the 
future. 
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